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Abstract. A time-dependent Kohn-Sliam scheme for ID bosons with contact 
interaction is derived based on a model of spinor formions. This model is 
specifically designed for the study of the strong interaction regime close to the 
Tonks gas. It allows us to treat the transition from the strongly-interacting 
Tonks-Girardeau to the weakly-interacting quasicondensate regime and provides 
an intuitive picture of the extent of fermionization in the system. An adiabatic 
local-density approximation is devised for the study of time-dependent processes. 
This scheme is shown to yield not only accurate ground-state properties but also 
overall features of the elementary excitation spectrum, which is described exactly 
in the Tonks-gas limit. 
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1. Introduction 

Density-functional theory (DFT) provides a unique framework for the treatment 
of quantum many-body systems beyond the realm of perturbation theory. While 
ground-state DFT, based on the theorems of Hohenbcrg, Kohn, and Sham [1,2] deals 
with the energy and the density profile of inhomogeneous ground states, the time- 
dependent generalization of these theorems by Runge and Gross [3] in principle allows 
the study of time-dependent processes and excited states. Further development of 
the time-dependent DFT has shown that it is advantageous to consider functionals 
simultaneously of the density and the current [4] in order to go beyond the simplistic 
adiabatic approximation and the microscopic Navier-Stokes equations for the electron 
fluid have been derived elegantly this way [o]. 

In this paper we will apply DFT to the study of fermionization in the system of 
ID bosons with contact interaction. In the case of infinite interaction strength, the 
model is also known as the Tonks- Girardeau gas [6] and maps exactly to the system 
of non-interacting spinless fermions with identical energy spectrum and single-particle 
density. In fact, it has been shown by Cheon and Shigehara [7] that the exact mapping 
between the original Bose system and a Fermi system can be extended to arbitrary 
interaction strength at the expense of a highly singular interaction in the fcrmionic 
picture. For the homogeneous Bose gas, exact solutions for stationary states at 
arbitrary interaction strength can be found using the Bcthe ansatz [(S,9, 10]. Although 
the exact solutions show a continuous transition between the perturbative regime 
of weakly-interacting bosons and the strongly correlated, fermionizcd regime, they 
do not prove very useful for the study of time-dependent processes or inhomogeneous 
situations. Neither do they provide us with a simple, intuitive picture of how strong the 
degree of fermionization is in a given system. The crossover from the quasicondcnsate 
to the fermionizcd Tonks gas has also been studied in Refs. [11, Vl] 

Inspired by the works of Haldane [13], Sutherland [14], and others on the concept 
of exclusion statistics describing a crossover between fcrmionic and bosonic statistics 
in ID systems, we study a model which explicitely allows this transition and provides 
an intuitive picture of the degree of fermionization. In order to devise a practical 
scheme for treating time-dependent and inhomogeneous systems we employ DFT and 
develop a time-dependent Kohn-Sham formalism based on the auxiliary model system 
of N non-interacting spin-(i^ — l)/2 fermions. This model is chosen because the spin 
degeneracy may simulate the level attraction or the bunching of single-particle quasi 
momenta in the interacting bosonic system [n]. The interaction energy of the Bose 
system is simulated by the kinetic energy of the spin-degenerate fermions. In this 
work we will study the model in the simplest and most generic approximation: the 
adiabatic local density approximation ( ALDA) . The limiting case of infinite interaction 
strength is obtained easily and is treated exactly with v = 1. The opposite limit 
of weak interaction can also be treated accurately within the proposed formalism 
with V = N where the perturbative Gross-Pitaevskii and Bogoliubov equations are 
recovered asymptotically. In the general case of arbitrary interaction strength, the spin 
degeneracy v is fixed by requiring the correct low-energy asymptotics of the excitation 
spectrum, which is analyzed in the framework of linear-response theory The resulting 
model is suited for the study of time-dependent processes in inhomogeneous ID Bose 
gases close to the Tonks-gas limit. The properties of this approximate model are 
analyzed and systematic improvements are suggested. 

Before developing the general theory we want to point out that DFT is not 
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commonly used to discuss the theory of ultra-cold bosonic gases although many 
standard approximations can be derived in this framework. Only few papers therefore 
explicitely suggested to apply DFT in this context, see e.g. Refs. [15,10,17,18,19,20,21], 
where the time-dependent DFT of superfluids [20,21] is particularly far developed. The 
case of strongly interacting bosons in ID has been considered by Kolomeisky et al. [17] 
who suggested a generalization of the Gross-Pitaevskii equation of the following form: 

in—ip{x,t) = I --^-Q^ +Vc^tix,t) +(j>{n{x,t)) ^ ip{x,t) (1) 

Here, 'ip{x,t) is a time-dependent complex field and n{x,t) = |'!/;(a;,i)|^ is the one- 
dimensional density. The function 4>{n) provides the nonlinear term. Kolomeisky et 
al. suggested to use the chemical potential of the Tonks- Girardeau gas 

<^(n) = 0(™)H^^n^ (2) 
2 m 

whereas the Gross-Pitaevskii equation is recovered for 

0(7i) = 0(«^)(n) = .giD,i. (3) 

We will refer to Eq. (1) as the bosonic LDA because this equation may be derived 
as a Kolm-Sham equation in the ALDA using a Bose condensate as the auxiliary non- 
interacting system for the Tonks gas and the weakly interacting Bose gas, respectively. 
The same equation also applies for ID bosons with arbitrary interaction strengths 
and was used by Ohberg and Santos [IS]. In this case, 4){n) represents the exactly 
known chemical potential of the Lieb-Liniger model (at T = 0). The bosonic LDA 
equation (1) is almost identical to a set of hydrodynamic equations [22, 2.'!] derived 
from general hydrodynamic arguments assuming local equilibrium, the only difference 
being a kinetic energy pressure term (see e.g. [24]). 

It may be questioned whether the choice of a Bose condensate as the auxiliary 
non-interacting system of the Kohn-Sham formalism is the ideal starting point for 
treating a strongly interacting Bose gas in ID which is not condensed [2''),(i]. Indeed 
it has been shown in the Tonks- Girardeau limit of impenetrable point bosons, that 
phase coherence is grossly overestimated by Eq. (1) leading to wrong predictions 
for interference properties [26,27]. Phase and density are conjugate variables in a 
quantum field theory and may coexist in the weakly-interacting Bose condensate, 
which is close to a classical field. It is known already in the perturbative regime, that 
phase fluctuations are very important in ID [28,12]. In the ultimate limit of the Tonks- 
Girardeau gas, the concept of phase looses its meaning and the governing equations 
should be rather based on the density alone. For this reason the hydrodynamic 
approximations employed in Refs. [22,23], where the concept of a phase docs not 
appear, are on a safer footing than the bosonic LDA equation of Kolomeisky where 
phase-coherence is built into the theory. 

Moreover, the excitation spectrum and density of the the Tonks- Girardeau gas 
is identical to the corresponding properties of a gas of independent spinless fermions 
due to the Bose- Fermi mapping. Already in his original paper [8], Lieb pointed out, 
that the fermionic character of the excitation spectrum prevails at finite interaction 
strength and introduced two elementary branches of excitations which he called type I 
and type II. Type I excitations result from exciting a particle from the edge of the ID 
Fermi sphere to an unoccupied orbital. Thus type I excitations are the particle-hole 
excitations of the highest possible energy for a given momentum. Type II excitations, 
on the contrary, are the excitations of the lowest possible energy for a given momentum. 



A density-functional approach to fermionization in the ID Bose gas 



4 



They can be seen as taking a particle from the inside of the Fermi sphere and putting 
it to the first free orbital on the outside. 

Only type I excitations are found by standard bosonic methods like Bogoliubov 
perturbation theory or linear response in Eq. (1) whereas type II excitations and the 
fermionic character of the spectrum remain hidden. We want to mention at this point 
that a possible connection between type II excitations and dark solitons in variations 
of Eq. (1) has been pointed out in the literature [29, 17,30,31]. 



2. An exotic Kohn-Sham scheme 



2.1. Variational principles 

The basic scheme of DPT is the Hohenberg-Kohn variational principle. It allows the 
ground-state energy and one-particle density n{x) of a given many-particle system 
to be found by variationally minimizing the energy functional 

E^[n] = Fyo^Iu] + j Vc^t{x)n{x)dx. (4) 

Here i'oxt(2;) is the external potential and the Hohenberg-Kohn functional i^HKi"-] is 
universal in that it does not depend on the external potential. The scheme of Kohn 
and Sham proceeds by partitioning the unknown functional i^HK [n] as 

FnK[n]=T,[n]+EKs[n], (5) 

where Ts[n] is the kinetic energy of a ficticious, auxiliary system of particles interacting 
only with a single-particle potential w^, which has the same density n(x) as the 
original system. The external potential Uext is contained in Vg — Vcxt + ''^ks where 
■fKsl"-] ~ ^E}^s/Sn is a universal functional of the density depending only on the 
structure and internal interactions of the original and auxiliary many-body problem 
but not on the external potential. Once an (approximate) expression for wksM i^^-s 
been found, the density of the system in a given external potential fcxt is easily 
found by solving the Kohn-Sham equations, a set of nonlinear independent-particle 
equations. Kohn and Sham in their original publication [2] used the auxiliary system 
of noninteracting spin-1/2 fermions while interested in systems of interacting electrons. 
However, it has been pointed out that the choice of the auxiliary system is somewhat 
arbitrary and rather a matter of convenience than dictated by laws of nature [32]. 



2.2. Hamiltonian 

We start with the usual Hamiltonian for N bosons in ID with point interactions of 
Ref. [fs], augmented with a possibly time-dependent external potential Vcxtix,t): 

( h"^ d 1 

i—l ^ ^ ^ i<j 

where giD = 2h^ / (mjaiDj) is the interaction parameter and any is the one-dimensional 
scattering length [11]. 

The exact many-body eigenstates and the complete excitation spectrum of the 
Hamiltonian (6) in the homogeneous system where Voxt = with periodic boundary 
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conditions have been found by Lieb and Liniger using the Bethc ansatz [6,;>]. The 
ground state energy is 

E'r'-N^eij) (7) 

where e{j) is the dimensionless energy per particle in the Licb-Liniger model and 
7 = 2/(n|aiD|) is the single dimensionless parameter of the homogeneous gas with the 
one-dimensional single-particle density n = N/L of N particles in a box of length L. 
The function 6(7) is defined as the solution of a Fredholm equation and can be obtained 
numerically to any desired accuracy. The chemical potential of the homogeneous Lieb- 
Liniger gas at density n is given by 

MLL(n) = -^E'r' - ^^/(7), (8) 
ai\ m|aiD| 

where fij) = [3e(7) — e'(7)]/7^ is a dimensionless function. The asymptotic behaviour 
of 6(7) and 7(7) is known for large and small 7 [f^] and, furthermore, the functions 
are tabulated for intermediate values of 7 in Ref. [22] |. 



2.3. Functionals and Kohn-Sham equations 

Consider a physical system of N noninteracting and indistinguishable fermions with 
spin {v — l)/2. For simplicity let us consider the case where v divides N and assume 
that all orbitals are nondegenerate. In the ground state N orbitals will be occupied 
with u particles each when vN = N . The kinetic energy is 

The Kohn-Shani equations arc found by functional differentiation with respect to 

Siifii (^) = ] ^ 'q^ + "cxt {x) + WKS [n] {x) > ipi (x) . (9) 

The potential wks [n] ~ ^^aw ^"^ mean-field potential of Kohn-Sham theory and 

represents the unknown rest of the functional. In regular electronic DFT, i;ks["] would 
be written as a sum of the Hartree potential and the exchange-correlation potential. 
The density is found self-consistently by summing over all occupied (g) orbitals 



n{x,t) = v^\ipi{x,t)\'^ 



The ground state density n{x) is obviously time-independent and only the N levels 
with the lowest energies will be occupied. The generalization to the corresponding 
time-dependent Kohn-Sham equations along the lines of the Runge-Gross theorems [?>] 
is obvious: 

ihdt^pi[x,t) = <^ ^^^^ +Vc^t{x,t) + VKs[n][x,t) > ifi{x,t). (10) 

\ The dimensionless chemical potential /{7) is known to have the expansion [s] 7(7) = 7r^7^^ — 
i^7-3 -I- 0(7"*) for 7 ^ 00 and 7(7) = 27-^ - §7"^ -I- 0(7^) for 7^0. We give the following 
convenient rational approximation /(7) « [7r^(48 + 27^ -|- 37^)]/{37[87r^ -1-8^^^-1-7(2 + 7)'^]}, which 
deviates with a maximum relative error of 0.12 from the numerical result at 7 2. It has the same 
expansion as given above for 7 — > 00. For small 7 it introduces an error of order 0(7^/'^). 
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We note that Eq. (10) is still exact. We now proceed in approximating WKsb^] ^ 
local-density approximation. 

To this end, the remaining unknown term i^xsM approximated by -EksM ~ 
E^^{n) = - T]'°"\ where T]^°"^ is the noninteracting kinetic energy of the 

homogeneous gas. It is simply found by summing the contributions from the orbitals 
in the Fermi sphere and multiplying with the degeneracy v. 

2m ^ \ L ) 

fe=-(W-l)/2 



where L is the length of a box with periodic boundary conditions, n ~ N/L = n is 
the line density and N has been assumed odd. 

Before taking the derivative '^"^g^ ^"^ in order to find an expression for vks we 
have to consider what to do with the term 1/7V. If TV = 1 the contribution from the 
kinetic energy vanishes and we arrive at the bosonic LDA of Eq. (1). However, if N is 
of the order of N, we may disregard the term 1/TV in the limit of large N. In fact it 
is not consistent to keep this term if we approximate the Lieb-Liniger energy with the 
bulk value as this is only valid in the limit N oo. In the following we will therfore 
drop this term and arrive at 

""KS^i^^t) = ^ (£^0°'" - ^kin) \n=n{x,t) 



where the chemical potential fihhin) of the homogeneous (Lieb-Liniger) gas at density 
n is given by Eq. (8). The remaining parameter v will be determined by a consistency 
condition on the excitation spectrum that will be derived in the following paragraph. 



2.4- Excitation spectrum from linear-response theory 

Lieb pointed out that the excitation spectrum of the ID Bose gas resembles that 
of the ID Fermi gas and likewise has a particle-hole type of structure — at any 
value of the interaction strength [9] ! The spectrum of the homogeneous gas may be 
equally well described by two branches of elementary excitations e^(fc) and e^^(fc). 
Both branches have the same slope for fc ^ and therefore lead to the same speed 
of sound Vs = de^/dk\k=o = de^^ /dk\k=o- Bogoliubov perturbation theory yields an 
approximation only for the branch e^(fc). In the time-dependent Kohn-Sham model 
introduced above, a particle-hole type excitation spectrum arises naturally. The 
elementary excitations may be found by solving for the resonant frequencies of the 
time-dependent equation (10) in the linear-response limit. The resulting excitation 
energies are given by differences of the orbital energies of the stationary Kohn-Sham 
orbitals modified by interaction contributions from the nonlinear mean- field Uxg'^. 

The derivation of the linear response equations is the same as for the LDA in the 
much studied electronic systems and very similar to the time-dependent Hartree-Fock 
or random-phase approximation (see e.g. [33]) and is summarized in the appendix. 
Here we want to note that the linear response correction to the spectrum of the Kohn- 
Sham system is proportional to the derivative dvKs/dn. From the explicit linear 
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response equations for the homogeneous system it can be seen that type I and type II 
excitations have a different slope at fc = 0, except for the situation where 

VKs(n) = 0. (13) 

This argument is detailed in the Appendix where the linear response equations in the 
small momentum limit are solved exactly. Using the expression (12) for vks we easily 
find the correct value for the number of modes in the fermionic Kohn-Sham equations 
(10): 

^ = ^\/-^/'(7)- (14) 

The value ofiV/A^^l/i^isa measure of the degree of fermionization in the system. 
When N/N ~ 1 the system is completely fermionized whereas N/N ^ indicates a 
purely bosonic system with a Bose condensate as the ground state. The functional 
dependence of N /N on the interaction strength 7 is shown in Figure 1. 

Qualitatively the picture agrees with the results of Girardeau and Wright [34] 
who modelled the BEC-Tonks crossover by a mixture of a pair-correlated Bose liquid 
with a non-interacting Fermi gas in a variational framework. We would like to point 
out, that iV is a global parameter whereas 7 is local and varies in an inhomogeneous 
system. Therefore the condition (14) should be used with care and we expect the 
present model to be most useful for nearly homogeneous systems. 

Although it would be desirable to find an approximation for the condensate 
fraction from our model, it should be noted that Kohn-Sham DFT as used here 
generally does not model the many-body wavefunction but gives direct access only to 
energies and the diagonal part of the single-particle density matrix. The information 
about off-diagonal long range order and condensate fraction, however, is contained in 
the inaccessible off-diagonal part of the single-particle density matrix. Nevertheless 
it should be pointed out that the fractional occupation of the Kohn-Sham orbitals 
A'^"^ in our model, which is reminiscent of a condensate fraction, vanishes in the 
thermodynamic limit where N —^ 00 ed, ^ = const, consistent with the general 
result that there is no Bose condensate in the T = homogeneous interacting Bose 
gas [35, 3G]. 

An important question related to the condensate fraction is local phase coherence 
[26,27]. From the construction of the Kohn-Sham model it appears difficult to define 
expectation values for the phase operator in this scheme. However, the time-dependent 
DFT does contain information about local phase coherence to the degree that is 
necessary to correctly describe the supression of interference patterns in the dynamics. 

The approximation scheme defined by Eqs. (10), (12), and (14) is the main subject 
of this paper and shall be referred to as fermionic LDA. In the remainder of this paper 
we will analyze the predictions of this model and compare with the exact solutions of 
the Lieb-Liniger model and other approximations. 

3. Properties 

3.1. Excitations of the homogeneous gas 

It is interesting to calculate the dispersion for the elementary excitations of type I and 
II, which can be done analytically since the linear-response corrections vanish due to 
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Figure 1. The parameter N/N from Eq. (14) indicates the degree of 
fermionization in the system. Here it is plotted as a function of the dimensionless 
interaction strength of the Lieb-Liniger model 7. 



Eq. (13). In the homogeneous gas, the mean field vks is constant and the Kohn-Sham 
single-particle energies become 



2t2 



-P 



where p is any integer and hkp 



hp2T:/L is the momentum. In the ground state 
all orbitals with p < F arc occupied where F — {N — l)/2 is the index of the Fermi 
momentum hkp. The energies of particle-hole excitations of momentum hkq are given 
by the difference 

fc2 



2m 



2 47r2 



(15) 



Figure 2 illustrates the discrete excitation spectrum obtained for a finite number 
of particles in a box with periodic boundary conditions. Type I and II excitations are 
the upper and lower bounds of the elementary excitation spectrum, respectively. The 
fermionic DFT is seen to slightly overestimate the energies of type I excitations and 
underestimate the energies of type II excitations at large momentum while the correct 
asymptotes arc obtained for small momenta. 

From Eq. (15) we can construct explicit expressions for type I and type II 
excitations. Type I excitations are defined by exciting a particle with the Fermi 
momentum hkp into an unoccupied orbital with p ^ F + q whereas type II excitations 
take a particle from inside the Fermi sphere to the lowest unoccupied orbital with 
p = F + 1. In terms of the dimensionless momentum 



P - 



hkq 

hn 



of the excitations we obtain 



2nq 













+p 







(16) 
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Figure 2. Shown arc excitation energies of a homogeneous system with 25 
particles at 7 = 16.3. Dots indicate the particle-hole excitation in the Kohn- 
Sham model with N/N = 0.8 from Eq. (15). The full and dotted lines show the 
exact type I and II excitations for this finite system from the numerical solution 
of Yang's equation [9, . 37], respectively. 



Correspondingly we find for type II 



(17) 



Expressions for the elementary excitation spectrum of the interacting gas can be 
derived using equation (14). In the thermodynamic limit wc obtain 



\P) 



h 

2m 



(P) = p 



V-273/'(7)+p 
V-2j^f'h)-p 



(18) 
(19) 



These dispersion relations ought to be compared with the exact ones and with earlier 
approximations. First of all we note that both dispersions relations indeed lead to the 
same speed of sound 

(20) 



dp 



TO aiD 



which is identical to the exact result of the Lieb-Liniger model [9] . 

The type II excitations, as discussed by Lieb, have the character of elementary 
excitations for p < n and at this point the slope de^^/dp vanishes. This condition is 
strictly fulfilled by Eq. (19) only in the Tonks-gas limit 7 — > 00 where the equations 
become exact. In the general case, e'^(p) of Eq. (19) has a maximum at Pm < tt and 
Pm approaches as 7 — > 0. This observation is related to the fact that the umklapp 
transition of taking one particle from one side of the Fermi-sphere to the other is 
not described correctly in the considered fermionic Kohn-Sham model. Therefore our 
approach is limited to describing type II excitations in the long-wavelength limit and 
works best for large 7. 
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3.2. Comparison with bosonic LDA 

The bosonic LDA of Eq. (1) is simply the single- mode case with TV = 1 of the general 
fermionic Kohn-Sham formalism in this paper. The linear-response equations (A. 2) for 
this case coincide with the Bogoliubov-de Gennes equations derived in Refs. [18,19]. 
It is instructive to derive explicit expressions for the excitations of a homogeneous gas, 
which appear not to be available in the literature. Eqs. (A. 2) now become 



(e<;-eoH —)Uq + —Vq = Tiu;Uq (21) 



N(f>\ v(f)' 

where (j)' = dcjy/dn and we have set X^p = Up-h and Y^p ~ Vh-p- With 
we find the simple solution 



Only the plus sign contributes here, the minus sign is a well understood artifact 
of linear response theory. In terms of the dimensionless momentum p ~ k/n we 
find the following final result for the excitation spectrum of the bosonic LDA with 
0(n) = ^i-Lhin)■. 

hcu = |lnVp2(p2- 273/'(7)). (24) 

In contrast to the excitation spectrum described by Lieb, there is only one excitation 
branch in the linear response of the bosonic LDA. Comparing with Eqs. (21) and (22) 
we see that the first terms in the expansions around fc = and k = oo are identical 
but there are discrepancies in between. The speed of sound is again the exact result 
Cs as a consequence of the compressibility sum rule which is obeyed by the LDA. 

We can recover the well-known Bogoliubov approximation in the limit of small 7 
by expanding 7(7) = 27"^ -I- 0(7"^/^) to find the usual dispersion relation 



{hio)Bos = —n^VPW+W)- (25) 
The speed of sound becomes 

Figure 3 shows the dispersion relations derived in this paragraph for a particular 
value of the interaction parameter 7. It can be seen clearly that the Bogoliubov 
approximation for the speed of sound is wrong, which is no surprise as 7 is not small 
in this example. However, it is also interesting to see that the dispersion derived from 
the bosonic LDA deviates significantly from the fermionic LDA. 
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Figure 3. Shown are excitation energies in the thermodynamic hmit at 7 = 16.3. 
The full and dotted lines are the excitations within the current Kohn-Sham 
approach of type I and II from Eqs. (18) and (19), respectively. The dashed 
line is the result from the linear response of the bosonic LDA (24). The dash- 
dotted line is the Bogoliubov dispersion (25) which shows the wrong slope at the 
origin because we are out of the perturbativc regime. 



3.3. Limits oj strong and weak interaction 

In the strongly-interacting limit we find v ~ 1 and N = N. In this case vks = and 
the Kohn-Sham equations become the equations for non-interacting spinless fermions. 
Due to the Bose- Fermi mapping theorem of Girardeau this gives the correct description 
of both the full excitation spectrum as well as the (time-dependent) diagonal single- 
particle densities. In particular we obtain the following equations for type I and II 
excitations: 

e\p) = —n^p[2TT+p] (26) 



,11/^ 



e 



{p)^—n'p[2Tr-p] (27) 
Zm 

It is interesting to note that the single-mode equation (1), which here is exactly the 
equation studied by Kolomcisky [17], does not have the correct limit of Eq. (26). 
Instead we find 



hio = nVp^(P^-27r2). (28) 

2m 

This limit underestimates the type I excitation energy. It has been suggested in 
Ref. [17] to link dark solitary waves in the single- mode equation to the type II 
excitations. The dispersion relation, however, gives only qualitative agreement and 
underestimates the energy at large momenta. We want to stress that our fermionic 
LDA, on the contrary, yields the exact excitation spectrum in the Tonks- Girardeau 
limit. 

In the limit of weak interactions we should obtain the single-mode equation with 
0(n) =A'll(7i), (29) 
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as was discussed in the eontext of the derivation of the Kohn-Sham equations. It was 
also discussed that the LDA approach for finite interactions followed in the main part 
of the paper cannot describe the limit of the Bose condensate with a single occupied 
mode correctly. Indeed it can be seen that Eq. (18) does not reduce to the Bogoliubov 
form (25) for small 7, although it has the correct low- and high-energy asymptotics. 
In order to obtain this limit correctly, a finite size correction should be added to the 
density functional. 

5.^. Static density profiles in the Thomas- Fermi approximation 

A consistency check for the static density profiles found from the fermionic LDA can be 
derived from applying a Thomas-Fermi approximation. The idea is to treat the Kohn- 
Sham model system of Eq. (9) as non-interacting Fermions (equivalent to a Tonks gas) 
in the external potential Uoxt {x) + {x) . Approximating now the density by the 
Thomas-Fermi approximation or hydrodynamic LDA as in Rcf. [22] leads to exactly 
the same equations for the density as the Thomas-Fermi approximation for the original 
system of interacting bosons discussed in Ref. [22]. This result is independent of the 
specific choice of N and provides a check for the consistency of the LDA between the 
different model systems. 

4. Conclusions 

In this paper we have devised a variational scheme based on Kohn-Sham DFT to 
calculate ground-state properties and time-dependent processes in the system of 
ID bosons with short-range interaction. The scheme will be particularly useful 
for the regime of strong interaction where fermionization is important and bosonic 
perturbation theory ceases to be applicable. In contrast to the previously suggested 
nonlinear Schrodinger equation of Kolomeisky [17], our theory has the correct strong- 
interaction limit for the density and for ground-state and excitation energies. The 
degree of fermionization in the model is determined by requiring consistency of 
the low-energy, low-momentum excitation spectrum. Employing the local-density 
approximation yields a parameter-free model. Interesting applications and test cases 
may arise in situation where coherence properties and long-wavelength excitations are 
important. Due to the construction of our model, it is best suited for situations with 
a large number of bosons and almost homogeneous densities. A specific example may 
be the study of shock waves in the homogeneous gas as studied by Damski [27, :^8]. In 
the decay of a shock wave, the usual bosonic LDA fails due to the fact that coherence 
is overestimated. Exact results are only available for the Tonks-Girardeau limit so 
far. For finite interaction strengths, our fermionic LDA can make predictions about 
the importance of coherences and numerical studies of this problem are under way. 
Another situation where the current theory can be applied is the ID Bose gas in a 
weak optical lattice. Interesting questions to study are the nature of excitations and 
the phase diagram as well as a comparison to the predictions of bosonic LDA. 

The DFT scheme studied in this paper is based on the careful choice of the 
kinetic energy functional determined by the Kohn-Sham auxiliary independent- 
particle system and the rest of the functional was approximated in the ALDA. The 
next logical step to improve performance of our model is to go beyond the ALDA and 
employ standard correction schemes based on current-density functional theory [4, ."i]. 
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Appendix A. Linear response equations 

Following the standard scheme of linear response theory [■'>'■'>] we introduce a time- 
dependent perturbation 

v,^t{x,t) = v,^t{x) + A[/+(x)e-'"* + /-(x)e-^'^*] 

in Eq. (10) where A is small. We look for solutions of Eq. (10) through first order in 
A. To zcroth order we obtain the stationary equation (9). Up to first order we expect 
to find solutions of the form 

Mx,t) = |(p.(x) + AX](^/.pe— * +y^y"*)^p(x)| (A.l) 

Here, the linear response of the occupied orbitals can be expanded in terms of 
unoccupied (^) orbitals because the time-evolution according to the Kohn-Sham 
equations preserves the orthonormality of the occupied orbitals to every order in A. 
Substituting the expression (A.l) into the time-dependent equation (10) and equating 
the first order terms with different time-dependence separately and finally setting 
= yields the following equations for small-amplitude free oscillations around the 
stationary state: 

i^hph'p'Xh'p' + Mhph'p'Yfi'p') = hujXhp (A. 2) 

h'e,p'^ 

i^lph'p'Yh'p' + Mlp^,,p,Xh'p') = -h^Yhp (A.3) 

In the basis of stationary Kohn-Sham orbitals [i.e. solutions of Eq. (9)] the matrices 
read 

Chph'p' = (ep - eh)Spp>6hh' + i^v'pf^,^,^ (A.4) 

Mhph'p' = vv'pp,,,,^, (A.5) 



where 



i]kl - / H'i V-^yf-jM 



(p*(x)(p*(x)^3^^ Lpk{x)ipi{x) dx. (A.6) 



dn 



i{x) 



A simple basis transformation may be applied to rewrite Eqs. (A.4) to (A.6) in a 
different basis, e.g. the coordinate representation. For the specific case of the bosonic 
LDA with iV = 1 we find the Bogoliubov-de Genncs equations of Rcf. [18]. In the 
limit 7 <C 1 the regular Bogoliubov-de Gennes equations of Gross-Pitaevskii theory 
are recovered. 

For studying the homogeneous system where v^-xt = and n is constant we 
introduce the usual box quantization with (j)ra{x) — e^^'^^ j \fL and kq = 2T:q/L. We 
find 

2m 
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and Vijki = 5i-kj-iv'^^/ L with uj^g = dv^s/dn. Due to translational symmetry 
the matrix equations block and it makes sense to label the blocks according to the 
momentum transfer hkp — hk^ ~ 2'n:qh/L. With defining -Cp p/ ~ ^{p-q)p-{p' -q)p' 
obtain 

^^p' ={^P~^P-,)6p..P' + '^ (A.7) 

Mly = ^ (A.8) 

Due to the presence of the Fermi sphere there are restrictions on the possible values 
of the indices as kp always has to be outside and kp-q has to be inside. 

We now consider the two cases q = 1, 2 explicitely in order to calculate the slopes 
of both the type 1 and type 11 excitation branches at zero momentum. For q = 1 
corresponding to the smallest possible momentum transfer, the matrices £ and M have 
one entry each and Eqs. (A. 2) and (A. 3) can be recast as a 2 by 2 matrix eigenvalue 
equation. This equation can be solved easily to yield the eigenvalues fiwi and —huji, 
where = + 2avv'Y^a,/L with a = N2tt^ / {mL'^i') and only the positive 

solution is physically relevant. Using Eqs. (12) and (8) we find that huji/pi — Cs 
where pi = 21:%/ L is the momentum transfer for this excitation and Cs is the speed of 
sound in the Lieb-Liniger model given by Eq. (20). Wc would like to point out that 
this result is independent of v and can be traced back to the compressibility sum rule 
which is granted by the LDA. 

We arc now going to check this result for the speed of sound by calculating the 
energies for the next larger momentum transfer with q ~ 2 where 2 particle-hole 
excitations give the first possibility to have different slopes for the type I and type II 
excitation branches. If both branches have the same slope and yield the same speed 
of sound, the energy difference between these two excitations should vanish in the 
thermodynamic limit. Wc have solved the eigenvalue equations (A. 2) and (A. 3) which 
are 4 by 4 with a specific symmetry. We omit the lengthy expressions but note the 
following observations: One of the solutions is independent of ly and gives the speed 
of sound as in Eq. (20) up to order 1/L. The other solution does depend on ly and for 
the difference between the squared slopes we find 

ih4'^r - it^^'^f] lv\ = £^<s + 0{\IL). (A.9) 

We see that this term remains finite in the thermodynamic limit unless uj^g vanishes, 
which gives a stringent criterion for the choice of = N jN as discussed in Section 2.4. 
If the criterion of vanishing w^g is violated, the model will not only develop phonon 
branches with different speed of sound but additionally mean-field instabilities may 
occur as indicated by complex eigenvalues of the linear response equations. In the 
thermodynamic limit of the homogeneous case this occurs whenever N is larger than 
the value determined by Eq. (14). 
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